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Incorporated with twisted boundary condition, Polyakov loop correlators can give a definition of 
the renormalized coupling. We employ this scheme for the step scaling method (with step size 
s = 2 ) in the search of conformal fixed point of SU(3) gauge theory with 12 massless flavors. 
Staggered fermion and plaquette gauge action are used in the lattice simulation with six different 
lattice sizes, L/a = 20, 16, 12, 10, 8 and 6. For the largest lattice size, L/a = 20, we used a large 
number of Graphics Processing Units (GPUs) and accumulated 3,000,000 trajectories in total. 
We found that the step scaling function a(u) is consistent with u in the low-energy region. This 
means the existence of conformal fixed point. Some details of our analysis and simulations will 
also be presented. 
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1. Introduction 

Gauge theories with many flavors are paid attention, as candidates of beyond the standard 
model physics which might give a solution to the hierarchy problem. In this work, we investi- 
gate a model which consists of SU(3) gauge theory coupled to 12 massless fermion flavors in the 
fundamental representation. We study the running of coupling constant in this system. Our main 
concern here is the existence of the infrared fixed point, which is a scheme-independent property 
of the model. In the infrared properties of such a system, nonpertubative effects are dominant, and 
the numerical simulations with the lattice regularization are very useful tools, as demonstrated in 
QCD. So far, this model has been studied by several groups using lattice simulations with different 
methods [[l]]-[f7|]. At this moment, they have not yet reached a clear consensus. In this work, we 
adopt the Twisted Polyakov Loop scheme [Q] |Q] [ JTOj ] to study the running coupling constant. This 
method is different from those used by others. We show result which supports the existence of 
conformal fixed point. 



We already published several articles using this method [Mil] - [ 14]. This work is the continua- 



tion of Refs. []12p[]13p. The data set used here is almost the same as our work reported in Ref. Ql lfl. 

This paper is organized as follows. In Sec. [2], we give a review on twisted boundary conditions 
and the Twisted Polyakov Loop scheme, which we used to calculate the renormalized coupling 
constant. In Sec. |3|, we explain the details of the set up of our simulations. The results are given in 
Sec. |j. Summary and comments are given in Sec.|5[ 

2. Twisted Polyakov Loop Scheme 

In this section, we explain the twisted boundary condition and our definition of renormalized 
coupling constant using Polyakov loops. The twisted boundary condition results in infrared cut-offs 
in both gluon and fermion propagators. Then, the calculation of step scaling function is feasible 
without using the Schroedinger Functional method, in which Dirichlet boundary conditions are 
required. 

The twisted boundary condition for the link variables is given as 

U fl (x+i>L v /a)=a v U li (x)al. (2.1) 

We apply this condition only for v = 1,2 direction. Here, = 1,2) are SU(3) matrices which 
satisfy the conditions D.1D.2 = e l2n ^ 3 D.2^i , = 1 » (^m) 3 = 1 an( ^ Tr l/V] = 0- Their explicit 

form can be [Q.i] a b = 8 a mo d(i+2,3)+i an d ^2 = diag(e~' 2jr / 3 ,e i27i: / 3 , 1). For v = 3,4 direction, we 
implement the periodic boundary condition, Q. 3 = Q. 4 = I. 
The Polyakov loops in the twisted directions are 



Pi(x 2 ,x 3 ,x 4 ) = Tr 



Y[ U i ( X l =j,X2,X 3 ,X 4 ) 
j 



The extra terms are added to maintain gauge invariance and translation invariance. The coupling 
constant can be defined as the ratio of the correlators of Polyakov loops in the twisted direction and 
periodic direction, 

(P 1 fa = 0) t fi(x 4 =W(2 a ))) _ 2 

(P 3 (x 4 = 0yP 3 (x 4 =L 4 /(2a))) 8 > { -> 
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with the factor k = 0.03184 . . . which is given by calculating one-gluon-exchange diagram. In the 
actual analysis, we average over data that are related by translation and 90-degree rotation, to reduce 
the statistical errors. To apply the twisted boundary condition to the system with fermions, we need 
to introduce another degree of freedom, "smell" , to ensure the single- valuedness of\y(x-\-L\\ + L2 2 ) 
under the application of two boundary twisting, in x\ and X2 directions, with the different ordering. 
The boundary condition is given by, 

V^(x + L v v) = eWa a v b x$ (n v )l a . (2.4) 

The indices a, /3 are labels for the "smell" degree of freedom. The factor e m ^ are multiplied to 
eliminate the zero-momentum mode in v = 1,2 direction. For v = 3,4 direction, we implement 
the periodic boundary condition, i/a(x + L v v) = \\f(x). The number of "smelF's is the same as 
the number of colors. The "smell" index does not combine with the gauge field, and it can be 
considered as the another flavor index. Under this condition, the number of flavors should be 
multiples of 4(taste) x 3(smell) for the dynamical simulations with staggered fermions without 
taking roots of fermion determinant. 



3. Setting of Simulations 

To study the scale dependence of the coupling constant by lattice simulations, we measure the 
step scaling function, which is defined as 



a(u,s) = g 



(3.1) 

«=* 2 tt) 



This is calculated by taking the continuum limit of the lattice step scaling function L(u,s,L/a). 

L(u,s,L/a) = gi ttice {p,sL/a)\ u=g 2^ il}L/a) , (3.2) 
a(u,s) = lim L(u,s,L/a) withLfixed. (3.3) 

«/L->0 

In this work, we fix s = 2, and define a(u) = a(u,2) and L(u,L/a) = £(w,2,L/a). 

We have done the hybrid Monte Carlo (HMC) simulation for SU(3) gauge theory with 12 
dynamical flavors by using staggered fermion. Twisted boundary conditions are applied for the 
gauge field and fermion field as explained in the previous section. For the gluon action, we 
use plaquette gauge action at /3 = 4.5 ~ 100. In all simulations, we set the mass of fermions 
to be zero. We used hyper-cubic box L\/a = Z^/a = Lj/a = L4/01 = L/a with the lattice size 
L/a between 6 and 20. Then the lattice step scaling function is calculated by using the com- 
binations of (L/a,sL/a) =(6,12), (8,16) and (10,20). The simulation was done by selecting the 
starting configuration to be in the true vacuum, which leads to a non-zero phase of Polyakov loop ( 
arctan(Im(P^)/Re(/^)) ~ ±27r/3 ) in ju = 3,4 direction, as discussed in Ref. [[I3|]. The simulation 
parameters are summarized in Table |l|. 

For the largest lattice L/a = 20, we perform the HMC simulations on GPUs. The simulation 
code was developed by using CUDA [|I5|]. We used around 50 ~ 100 GPUs at the same time. 
Most of them are Tesla 1060. We used one GPU for one Markov chain. For each parameter set, 
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Lattice Size 


j3 


6 4 


4.5, 4.7, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0 , 50.0, 100.0 


8 4 


4.5, 4.7, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0, 50.0, 99.0 


10 4 


4.5, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 50.0, 99.0 


12 4 


4.5, 4.7, 5.0, 5.3, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0, 50.0, 99.0 


16 4 


5.3, 5.5, 5.7, 6.0, 6.5, 7.0, 7.5, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0, 50.0, 99.0 


20 4 


5.7, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0, 50.0 



Table 1: Summary of simulation parameters 



we produced 10 to 50 Markov chains at the same time, starting from the same seed configuration 
which is thermalized for that parameter set. To remove the correlation between the different Markov 
chains with the same starting configuration, we discarded 500 trajectories in the beginning of each 
Markov chain. The performance of the simulation with GPU is around 25 GFlops on average 
(sustained). 

To extract o{u,s) in Eq(3J), in practice, we have done simulations with all beta values listed 
in Table [j] first, then we perform the interpolation in j3 for the renormalized coupling at fixed L/a. 
The lattice step scaling functions (3.2) are calculated from the interpolated values of the renormal- 
ized coupling constant. For this interpolation, we used non-decreasing polynomials with degree 
(2?i + l), 



f{x)= J dx(^c,x'\ =h + h 1 x^ h h 2 n+ix 2n+l . 



(3.4) 



with x = 1/(5. Here, we choose ho to be and h\{= c 2 ,) to be 6 to match with the perturbative 
expansion. Thus the number of free parameters becomes (n— 1). 



4. Results 



In Fig. [I], we show the /3 -interpolation of renormalized coupling constant. The degree of the 
polynomial (2n+ 1) in Eq. ( |3.4| ) is chosen, at fixed L/a, such that it leads to the smallest % 2 / (d.o.f.), 
as shown in Table^. From Fig. [|-a to Fig.[I]-c, it is clear that the smaller L/a is, the larger the 
difference between g 2 Rtt (fi ,2L/ a) and g 2 att (p,L/a) is. This is particularly obvious in the infrared 
regime gf att > 2. It means that, for a given input coupling u, the lattice step scaling function is 
smaller for the larger L/a, as shown in Fig. ^-a. From Fig.[j]-d, we see that gf attice (p,L' /a) > 
lattice (PiL/a) for L' /a > L/a except for ultraviolet region. This means, outside the asymptotic- 
freedom regime, the lattice step scaling function ~L{u,s,L/a) is larger than u for any step size s. 

In Fig. ||, the lattice step scaling function and the continuum limit are shown. From the Fig. ||- 
b, we see that linear fit gives much smaller # 2 /(d.o.f) in the region u = 1.8 ~ 2.4, on the other 
hand, constant fit gives a larger value of # 2 /(d.o.f) ~ 0.8 to 20 in that region. For u = 1 ~ 1.7, 
the % 2 / (d.o.f) from constant fit and linear fit are about the same. From these, we conclude that the 
linear fit is much preferable than the constant fit. In Fig.^|-c, we show the growth ratio o{u)/u. 
To estimate the systematic error, we changed the number of fit parameters in the j8 interpolation 
Eq. (3.4), from the choices given in TableBby ±1 for all lattice sizes. We calculated the step scaling 
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Size L/a 


6 


8 


10 


12 


16 


20 


Number of points 


17 


17 


19 


18 


17 


13 


Number of fitting parameters 


8 


8 


6 


5 


8 


7 


Degree of polynomial 


19 


19 


15 


13 


19 


17 


Z 2 /d.o.f. 


1.04 


0.97 


0.59 


1.65 


0.54 


0.61 



Table 2: j 2 /d.o.f. for the fitting of beta. 



function using all possible 3 6 combinations for this interpolation. The narrow (red) band in Fig. ||-c 
is the statistical error associated with the procedure using the best fit (with parameters in Table |f) 
for the /3 interpolation. On the other hand, the wide (purple) band indicates the largest deviation 
(including the statistical fluctuations) from the central value of the best fit in all the j3 interpolations 
that we carried out. Therefore, this band contains both statistical and systematic errors. In Fig.^|-c, 
it shows that o{u)/u is consistent with one in 1.79 < u. This supports that this theory contains an 
infrared fixed point. 



5. Summary and Comments 

In this work, we show the existence of the infrared fixed point (IRFP) in SU(3) gauge theory 
with 12 massless fermions in the fundamental representation, by using the step scaling method with 
the Twisted Polyakov Loop scheme. The results presented in this article are obtained with step size 
s = 2. The location of IRFP reported here is consistent with our work [ |Tl| ] using step size s = 1.5. 

The lattice step scaling function exhibits L(u,L/a) > u at finite lattice spacing, except for 
the ultraviolet region. On the other hand, the step scaling function in the continuum limit shows 
that o{u) ~ u in the infrared region. This property at finite lattice spacing is different from the 
results in the step scaling function study using the Schroedinger Functional scheme [JTJ] [^]] , in which 
L(u,L/a) ~ u in the infrared region even at the finite lattice spacing. 
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Figure 1: (a, b, c) - The renormalized coupling from the simulation are shown by points with error bars. Fit 
functions (± error by A% 2 = 1) are shown as curves. The lattice sizes L/a are (6, 12), (8, 16) and (10,20), 
respectively. The value of input coupling u for the step scaling function at the left and right edges in Fig. 2-c, 
0.6 and 2.4, are shown for the guide of eyes, (d) - Central values of the fitting function for all sizes are 
shown. 
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Figure 2: (a) - The lattice step scaling function T.(u,L/a) and the step scaling function in the continuum 
limit a{u) with linear fit. (b) - The quality of fit % 2 /{d.o.f) for linear fit, constant fit and constant fit with 
only finer lattice, (c) - Growth ratio a{u)ju with linear fit. Systematic errors are given by changing the 
degree of polynomial in fit function. The result of perturbative calculation using one-loop and two-loop are 
also shown. — Step size is 2 for all panels (a),(b),(c). 
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